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We present results for lattice QCD with staggered fermions in the limit of infinite gauge coupling, 
obtained from a worm-type Monte Carlo algorithm on a discrete spatial lattice but with contin- 
uous Euclidean time. This is achieved by sending both the anisotropy parameter 7^ ~ 0/0^ and 
the number of time-slices Nr to infinity, keeping the ratio ~ aT fixed. In this limit, am- 

biguities arising from the anisotropy parameter 7 are eliminated and discretization errors usually 
introduced by a finite temporal lattice extent A^^ are absent. The obvious gain is that no contin- 
uum extrapolation Nr ^ °° has to be carried out. Moreover, the algorithm is faster and the sign 
problem disappears completely. 

As a first application, we determine the phase diagram as a function of temperature and real and 
imaginary baryon chemical potential. We compare our computations with those on lattices with 
discrete Euclidean time. Discretization errors due to finite A^^ in previous studies turn out to be 
large at low temperatures. 
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1. Introduction 

The determination of the QCD phase diagram, in particular the location of the critical point, is 
an important, long standing problem, requiring non-perturbative methods. In lattice QCD, several 
approaches have been developed to investigate the phase transition from the hadronic matter to the 
quark gluon plasma, but all of them are limited to Hb/T < I, with /ig the baryon chemical potential 
[1]. The reason for this is the notorious sign problem, which arises because the fermion determinant 
for finite /ig becomes complex, and importance sampling is no longer applicable. In QCD, the sign 
problem is severe. The relative fluctuations of the complex phase factor grow exponentially with 
the lattice volume. However, in the strong coupling limit of lattice QCD (SC-QCD) discussed 
below, the order of integration is reversed: the gauge links are integrated out first, and the partition 
function is expressed as a gas of hadron world lines. There is no fermion determinant, and the sign 
problem is much milder. This allows us to obtain the full (/Xg, T) phase diagram. 



2. Strong Coupling Lattice QCD in the Continuous Time Formulation 

In SC-QCD, the gauge coupling is sent to infinity and hence the coefficient of the plaquette 
term j8 = is sent to zero. Thus, the Yang Mills part F^yF^y of the action is absent. Then, 
the gauge fields in the covariant derivative can be integrated out analytically. However, as a con- 
sequence of the strong coupling limit, the lattice spacing a becomes very large, and no continuum 
limit can be considered. The degrees of freedom in SC-QCD live on a crystal. We study the SC 
Umit for one flavor of staggered fermions.^ The action is given by the fermionic part only: 
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(2.1) 



with niq the quark mass and /i = i/ig the quark chemical potential. The anisotropy pai^ameter 7 
modifies the Dirac coupling of the temporal part. It will be discussed in detail below. 

Following the procedure discussed in [2], the gauge link integration for gauge group SU(A'^c) 
can be performed analytically, as the integration factorizes in Eq. (2. 1): 

J X n kf,{x)=0 ^^c-K^{X). \\ y J 

+ K {p{x,yf^B{x)B{x + fL) + {-p{y,x)f^B{x + fl)B{x)) } (2.2) 
M{x) = x{x)x{x), B{x) = ^^ ei^,„i,^Xh W • • • Xi,,. W (2-3) 

The new degrees of freedom are mesons M{x) and baryons B{x): The first part of the expression 
(2.2) describes the mesonic sector (where k^{x) counts the number of meson hoppings), and the 



'in the continuum limit of staggered fermions, one flavor represents four degenerate tastes. In the strong coupling 
limit, one simply has one flavor of spinless fermions on a crystal. 
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second part describes the baryonic sector, which involves a /i-dependent weight p{x,y). After 
performing the Grassmann integrals analytically, the strong coupling partition function is: 

^K,M.)= I n ^^^j^r^'^'^'U^^^^ (2-4) 

with the constraint: nx+ ^ k^{x)=Nc, Vx G Vm (2.5) 

/i=±6,...±3 

This system, obtained by an exact rewriting with no approximation other than the strong coupling 
hmit involved, can be described by confined, colorless, discrete degrees of freedom: 

• Mesonic degrees of freedom: monomers «v € {0, . . .Nc} anddimers kjx{x) G {0, . . . A'f } which 
are non-oriented meson hoppings. They obey the Grassmann constraint Eq. (2.5) on mesonic 
sites Vm C N(j^ xNr- 

• Baryonic degrees of freedom: they form oriented, self-avoiding loops / with weight w(/) 
involving the chemical potential, the sign a{l) = ±1 and winding number r(/) which depend 
on the geometry of the loops. Baiyonic sites are not touched by mesons: VbOVm = Na^ x A^^- 

Here, we consider the chiral limit, niq = 0. Then, from Eq. (2.4), monomers are absent: rix = 0. 

In Eq. (2.1) we have introduced an anisotropy 7 in the Dirac couplings. This complication 
is necessary because the chiral restoration temperature is given by roughly ~ 1.5, and on an 
isotropic lattice with aT = l/Nt we could not address the physics of interest. Moreover, with the 
plaquette term being zero, varying 7 is the only way to vary the temperature continuously. The 
temperature, given by the inverse of the lattice extent in the temporal direction, is thus 

r = 4r with f{Y) = a/a,. (2.6) 

However, the functional dependence /{y) of the ratio of the spatial and temporal lattice spacings 
on 7 is not known. Naive inspection of the derivatives in Eq. (2.1) would indicate /{y) = 7, but this 
only holds at weak coupling. The mean field approximation of SC-QCD for SU(A'^c) gauge group 
based on 1 / J-expansion with d the spatial dimension [3] yields for the critical anisoti'opy 

^ ' 6(iV, + 3) ' ^ ' 

suggesting that aTc = is the sensible, A^^- -independent identification in leading order in l/d. 
This is also confirmed numerically: In Fig. 1 (top left) we show the variation with A^,- of the chiral 
susceptibility close to the U(3) chiral transition, and in Fig. 1 (top right) we show the vaiiation with 
Nt: of the U(3) transition temperature defined as iljNx- Both quantities approach continuous time 
(CT) limits, but not monotonically. Corrections to the CT limit, e.g. for the critical temperature Tc, 
can be parameterized as 

ar,(A^O = «7f ^ + fi/^T + C/^T^ + ^(A^r^) (2.8) 



^Note that baryons transform under gauge transformation e f/(3) as B[x) B(x) detii, hence they are not U(3) 
gauge invariant. U(3) describes a purely mesonic system ()f = in Eq. (2.2)), while SU(3) contains baryons (k: = 1). 
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where B and C have opposite signs. To circumvent this difficult extrapolation problem, we have 
developed an algorithm which samples directly the CT partition function in the limit 

A^^^oo, 7-^00, f/Nj = aT fixed. (2.9) 

Hence we are left with only one parameter setting the thermal properties, and all discretization 
errors introduced by a finite A^t are removed. Additionally, an algorithm operating in this limit 
has several advantages: There is no need to perform the continuum extrapolation A^^- — > 00, which 
allows to estimate critical temperatures more precisely, with a faster algorithm. And in the baryonic 
sector of the partition function great simplifications occur: Baryons become static in the CT limit, 
hence the sign problem is completely absent. 

We now explain the partition function we have used for Monte Carlo in the CT limit: after 
factorizing Eq. (2.4) into spatial and temporal parts, the spatial part simplifies greatly when taking 
into account which configurations are suppressed with powers of 7^^: 

The identity holds only approximately if 7 is finite, but becomes exact as 7 — )• 00 because spatial 
dimers with multiplicity fc,(x) > 1 are suppressed by powers of 7 and hence disappear when 7—)- 00. 
This is illustrated in Fig. 2 (right): as the temporal lattice spacing ~ c?/7^ — )■ 0, multiple spatial 
dimers become resolved into single dimers. The overall number of spatial dimers remains finite in 
the CT limit, as the sum over ^(7^) sites compensates the 1/7^ suppression. Temporal dimers can 
be arranged in chains of alternating ^-dimers and {Nc — ^) -dimers. In particular, for Nc = 3, we 
denote 3-0-chains as dashed lines, and 2-1 -chains as solid lines (see Fig. 2). 

The crucial observation is that the weight of these chains in the partition function is indepen- 
dent of their length, because the weight of each fc-dimer is the inverse of that of the {Nc — ^)-dimer. 
Hence, the weight of a configuration will only depend on the kind and number of vertices at which 
spatial hoppings are attached to solid/dashed lines, not on their positions. 

For SC-QCD with Nc = 3, there are two kinds of vertices, "L"-vertices of weight vl = 7^\ 
where dashed and solid lines join, and "T"-vertices of weight vj- = ly^^ /\/3, where a solid line 
emits/absorbs a spatial dimer. The partition function can now be written in terms of these vertices: 




VL = 1, VT = 2/V3. (2.11) 



The exponents ?ii(x) and nT{x) in Eq. (2.11)denote the number of L- vertices and T- vertices at spa- 
tial site X G Vm- In contrast to meson hoppings, spatial baryon hoppings are completely suppressed 
in the CT limit by a factor 7^/7^ = 7^- Hence, baryons are static in continuous time and the sign 
problem has completely vanished! In the limit N^ ^ °° the CT partition function can be written 
simply in terms of j8 = A^t/7^> up to some normalization constant: 

^CT{li)= l,{P/2f e'^^Pv^^ with k= h, NT = Z^Tix) (2.12) 

kem -^'en b=(x,i) ^ 
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Figure 1: Top left: continuous-time limit of the chiral 
susceptibility x at Y^/Nr = 1.8 in the purely mesonic 
system U(3), exhibiting typical non-monotonic be- 
havior in l/Nr. Note that the suppression of double 
and/or triple spatial dimers has no effect on the CT 
limit. 7b;? right: Non-monotonic continuum extrap- 
olation for aTc and comparison with the continuous 
time limit. Bottom: Reweighting of the energy based 
on the partition function Eq. (2.12). 
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where Ft is the set of all configurations W equivalent up to time shifts of the vertices and with a 
total number k of spatial hoppings, and B is the baryon number. In Fig. 1 (bottom) we illustrate 
the single histogram reweighting of the energy in j8 based on this partition function. An important 
property of the above partition function Eq. (2.12) is that spatial dimers are distributed uniformly 
in time. The lengths Aj8 of dashed or solid time intervals (which are related to the number of L- 
and T-vertices) are then, according to a Poisson process, exponentially distributed: 



P(A/3) ~ exp(-AAj8), A/3 G [0,j3 = l/aT] 
X=dM{x)/4, dM{x) =2d — ^nB{x±fi) 



(2.13) 
(2.14) 



with A the "decay constant" for spatial dimer emissions. The presence of baryons results in A 
being space dependent, with dM{x) the number of mesonic neighbors at a given coordinate x. Non- 
trivial meson correlations arise from the entropy of the various configurations. Likewise, baryonic 
interactions beyond the original hard core repulsion are due to the modification they induce on the 
meson bath, and thus also arise from entropy. 



3. Continuous Time Worm Algorithm 

Continuous time (CT) algorithms ai^e now widely used in quantum Monte Carlo (see e.g. 
[4, 5]), but to our knowledge have not yet been applied to quantum field theories. Special difficulties 
associated with the local gauge symmetry are absent in our case, since gauge fields have been 
analytically integrated out. The CT algorithm used here is a Worm-type algorithm, similar to 
the directed path algorithm introduced for SC-QCD in [6]. In the latter, sites are partitioned into 
"active" and "passive" depending on their parity with respect to the worm tail position. Here, we 
similarly have "emission" and "absorption" sites, giving the spatial dimers a consistent orientation. 
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Figure 2: Le/f; updating rules for the continuous 
time algorithm. Top: illustrative 2-dim. configurations 
(time flows to the right), in discrete time {leji) and con- 
tinuous time {right). Note how the latter lacks multi- 
ple spatial dimers and has only static baryons lines. 



In Fig. 2 we outline the updating rules of the continuous time worm, and also show a typical 2- 
dimensional configuration in continuous time. 

In contrast to simulations at finite N^, there is no need for a baryonic worm update. A simple 
heatbath update on static lines with no spatial dimers attached is sufficient: positively (negatively) 
oriented baryons are (dis)favored by a factor exp(ib3/i /T) over static meson lines, and their weight 
is always positive. We found it very useful to resum baryonic and static mesonic degrees of free- 
dom into so-called static polymers similarly to the MDP algorithm [7] : first, this allows to extend 
simulations to arbitrary imaginary chemical potential. Second, this enables us to adapt the Wang- 
Landau method for determining the nuclear transition at low temperature very accurately. This 
technique reduces the uncertainty on the first order line pLc{T) greatly. The polymer formulation we 
have used resums static mesons, baryons and anti-baryons. A configuration with P static polymers 
has weight^ w{P) = (4 + 2cosh(3/i/r))^. In our simulation, we only keep track of the polymer 
number P = Dq + B^ + B , which is the sum of the numbers Dq, , and B of static mesons, 
baryons and anti-baryons, respectively. We calculate the baryon number B = B^ —B from P via a 
trinomial distribution. In the Wang-Landau method, we hence obtain at fixed T the density of states 
for baryon number B, g{B, T), from the density of states for static polymer number P, g{P, T). 

4. Results on the SC-QCD Phase Diagram 

In SC-QCD at low temperature, chiral symmetry, i. e. the i7A(l) symmetry of the one-flavor 
staggered action, is spontaneously broken according to Ui{\) x Ur{\) — )■ Uv{\-) and becomes re- 
stored at some critical temperature Tc{pL). Our new results for the phase boundary Tc{pL) in the 
chiral limit lUq = eliminate systematic errors affecting previous findings based on mean field 
approximations [8] or Monte Carlo for fixed A^^ [9]. As previously found, the phase transition is 
second order at small /i, and first order at low T . The tricritical point is located close to its ear- 
lier Nr = 4- estimate [9], if one uses a/a^ = y^, i.e. aT = y^/Nr, ajx = y^{atlJ-)- However, the 
re-entrance predicted by mean field analysis [8] (which fixes 7=1 and varies A^^ G M) and seen in 
earlier Monte Carlo studies [9, 10] is absent, see Fig. 3 (top left). As illustrated in Fig. 3 (right) 
the discretization error from A^^ = 4 or 2 becomes very large at low T, and mimics re-entrance. 
We have also considered an imaginary chemical potential. As shown in Fig. 3 (bottom), the phase 

^Note that there are 4 kinds of static mesons: dashed hnes and solid hnes of even and odd parity. 
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Figure 3: Top left: Comparison of SC-QCD phase 
diagram obtained from mean field [8] and in contin- 
uous time. The former predicts a re-entrance which 
is not present in CT Monte Carlo. Top right: The 
SC-QCD phase diagram obtained with Nt = 2,4 us- 
ing aT = Y^/Nr, all = /^(a^/i) [9] and in continuous 
time. Bottom: Extension of the phase diagram in- 
cluding imaginary chemical potential, with no signal 
for a Roberge Weiss transition (x=-l) at high T. 



diagram displays Z(3) periodicity, but the imaginary /i transition remains second-order and there 
is no Roberge-Weiss transition at high temperature. This is natural: at high T the partition function 
Eq. (2.12) is dominated by the ^ = term which is analytic in jj.. 

In summary, the continuous time formalism allows for a final, unambiguous determination of 
the strong coupling phase diagram, and is suitable for further extensions like a second quark flavor. 
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